Semiclassical analysis of the Nonequilibrium Local Polaron 
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A resonant level strongly coupled to a local phonon under nonequilibrium conditions is investi- 
gated. The nonequilibrium Hartree-Fock approximation is shown to correspond to approximating 
the steady state density matrix by delta functions at field values to which the local dynamics relaxes 
in a semiclassical limit. If multiple solutions exist, all are shown to make nonvanishing contributions 
to physical quantities: multistability does not exist. Nonequilibrium effects are shown to produce 
decoherence, causing the standard expansions to converge and preventing the formation of a polaron 
feature in the spectral function. The formalism also applies to the nonequilibrium Kondo problem. 
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The quantum mechanics of nonequilibrium systems is 
a subject of fundamental importance and of great cur- 
rent interest, for example in the context of prospective 
'single molecule devices i In equilibrium problems, non- 
perturbative analysis based on solutions of Hartree-Fock 
equations (which may be understood as saddle points of 
functional integrals) has led to important insights. In 
nonequilibrium problems, mean field equations may be 
formulate d^i'^'^ and indeed exhibit a richer structure than 
in the corresponding equilibrium problems, but it has not 
been clear how to select the relevant solutions or to sys- 
tematically compute corrections. 

In this paper we investigate a simple model which in- 
dicates a resolution to these issues. We find that the 
relevant solutions are selected by the steady state den- 
sity matrix, which in a semiclassical, weakly nonequilib- 
rium limit is found to become very sharply peaked at 
field values corresponding to local minima of a 'pseu- 
doenergy' which we define. The formalism also shows 
that if multiple minima exist, all contribute, with weights 
varying smoothly as parameters change and shows how 
departures from equilibrium lead to a decoherence which 
suppresses characteristically quantal effects such as the 
formation of a polaron resonance. 

We consider a single level which may be occupied 
by or 1 spinless electrons (creation operator d^) and 
is coupled to two leads j — L, R (creation operators 
a^j,). The leads are assumed to be reservoirs specified 
by chemical potentials fij and inverse temperature (3j 
(which we typically set to (3j = oo). The electrons 
interact with an oscillator (dimensionless displacement 
coordinate q and momentum p) of mass Mph and en- 
ergy U{q)). This is a local version of the familiar 'po- 
laron problem'^- and captures an important aspect of the 
physics of prospective single- molecule devices^'^. The 
Hamiltonian is H ^ Hce + Hmn + Hdot + Hph + Hj; 
with 

Hce = ^ £k,aal jak,j; Hdot = eod^d; (1) 
Htun = ^ {d^akj + h.c) (2) 



Hi^XqSd- Hph^-^ + Uiq) (3) 

At Mph = oo the model is analytically solvable, with 
properties determined by the Green functions 

S',':,'0 w 

where g{t) = J |^e-''^*g(w) with 

= ZJ^T^^x^^^i^J+IS) ^ 9qAiujy (5) 
5gA'(t^) =2i^j.^^^^[aj(tj)tanh(/5j(w-^j))], (6) 

and S = E^+Sfi, ^j{z) = Y,p z-e\ = ImT,j{uj- 

iS)/ ({uj - £0 - A(j - ReJ:{Lj)f + /m£2(w)^ . It is also 
useful to define b — Regu and the 'pseudoenergy' 

m = U{q) + fdq'Xq' [ K - 9,^. + 9, A 
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which in equilibrium becomes the q-dependent part of 
the total energy. Fig. shows a possible form of $(9). 

The equilibrium physics of this model is very well un- 
derstood, and is conveniently viewed as a functional inte- 
gral over trajectories q{t). In the Hmits Mph,P 00 the 
integral is dominated by those paths for which q takes 
the definite value qa minimizing the energy. At finite 
Mpfi, other paths become important: two crucial classes 
are the 'gaussian fluctuation' paths involving small ex- 
cursions (characteristic frequency uJa) from the miniman 
and, if ^eqil) has the form shown in Fig ^ and the bar- 
rier height H is sufficiently large, tunnelling (instanton) 
processes during which q goes rapidly from the vicin- 
ity of the global minimum (Ei in the notations of Fig 
^ to the higher energy minimum (here E2), spends a 
time of order AE^^ = [E2 — Ei) near the higher min- 
imum, and then returns to the vicinity of the global min- 
imum. If H is sufficiently large, then tunnelling pro- 
cesses are rare, but when E2 ~ Ei, the time spent in 
the higher minimum becomes longer than the interval 
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FIG. 1; Assumed dependence of 'pseudoenergy' (Eq. on 
oscillator coordinate q. 

between tunnelling events and the dilute instanton ap- 
proximation breaks down, signalling the formation of a 
polaron resonance in the density of states. In this paper 
we investigate the changes occurring when the system 
is driven out of equilibrium by application of a chemi- 
cal potential difference A/x between the two leads. We 
present a general discussion but focus most attention on 
the nonequilibrium polaron limit H » uja » A/i; with 
Pj — oo but AE/Afi arbitrary. In our actual calcula- 
tions we also assume that departures from equilibrium 
are small enough that we may neglect density-of-states 
variations: {Afj,)dlna/duj < 1. 

To analyse the out of equilibrium behavior we use 
the standard^ extension of Feynman diagrammatics to 
the Keldysh contour which consists of a time-ordered 
(— ) branch extending from t = —oo to t = +oo, fol- 
lowed by an anti-time-ordered branch extending from 
-|-oo back to — oo. The Keldysh diagrammatics may be 
obtained by functional differentiation, with respect to 
source terms ria(t), of a generating functional W^[{?7a] 
which can be formulated^ as a coherent state path in- 
tegral: 

W[{Va}]= J dq+dq^p{q+,q^) J Vq(t)SK[-oo+, -oo^] 

(8) 

Here the functional integral / T>q{t) is over all paths be- 
ginning at (7_ at t = — oo on the time ordered contour and 
ending at q+ ait = — oo on the antitime ordered contour, 

while SK[tut2]Va{t)] = rKe'-f'^'**"(*)^[''"(*)l the time 
evolution operator on the Keldysh contour and a{t) = ±1 
according to whether the time is on — or + branch. The 
contributions of paths with endpoints (7+ , g_ are weighted 
by the appropriate element of the steady state density 
matrix p which is the long time limit of an equation which 
may be expressed in path integral language as 

p{q+,q-;t;Tl) = jdq'+dq'_jV'q'{t') (9) 

SK[t+, to+\p{q'^ ,q'_\to)SK [to- , ] 



Here the / T>'q'{t') is over all paths which begin at q'_ 
on the lower contour at time to and end at g_ time t 
on the lower contour, along with the time reversed paths 
which begin at q'_^_ (t) on the upper contour and return to 
g+ at time to- We require the long time behavior after 
transients have decayed. We do not find non-steady long 
time limits such as limit cycles or chaos: p in Eq|Slis the 
time independent solution of Eq|51 

Because eqs I1I2I3I involve a finite system coupled to 
two infinite reservoirs, the trace over electron opera- 
tors may be performed^. The combination of time evo- 
lution operators needed in eqs. I8I9I may be written 
in the interaction representation as TZ{t2,ti; {q{t)}) = 
{TKExp[-i dt{Xqc {n--n+) + Xqqg{n^+nj^))]) where 
q^lt) = i^M^kw^ q^^t) = i-to±m^ denotes expec- 
tation value in the reservoir defined by the leads (note 
this trace depends on the entire trajectory q{t)) and 
= d^d on the indicated contour. In the physical 
problem the quantum field coupling constant Xg = X. 
Differentiation with respect to Xq and application of the 
usual linked cluster arguments leads to 

7^(^2,tl;fe(t),g,(^)}) = e--fn (lo) 

Here Cph+s is the phonon action Mph{dtqqdtqc) — 
^j-9g+9c-| _ jj^ qg-Qc -| s^ppienjented by source terms, 

C{t)=qq{t)[ dXqTr[a^G^^{t,t)] (11) 
Jo 

and Gg> solves 

Ga, =g + g*v*GA, (12) 

with v(t) = {Xqc{t)l + Xqqq{t)ax)- We have verified that 
an expansion in powers of A about the weak coupling limit 
A = reproduces results obtained^ by standard Keldysh 
diagrammatic analysis. 

We now turn to the semiclassical analysis. In the large 
Mph limit, one expects stationary (time-independent) 
paths to dominate the physics. However, in the Keldysh 
formalism, stationary paths have = 0, so in the ab- 
sence of source terms the time evolution operator in Eq 
|Hlis unity for all stationary paths, providing no basis for 
selection. Instead, the important paths are selected by 
the density matrix; i.e. from the solution of Eq|3 

To understand the dynamics implied by Eq|^we first 
assume that p{q,q') is strongly peaked near q = q' = qa- 

We may then expand <& = £^0 + —^5—^ (<? — 9a )^ and anal- 
yse Eq 1^1 by the usual perturbative methods^'. In the 
large mass, weakly nonequilibrium limit wc find relax- 
ation, with a rate of order X'^uJa{aL -t- oj^)^ to a sum of 
functions sharply peaked about the values which min- 
imize the pseudoenergy. 

a 

To leading order in the electron-phonon coupling 
A we find r,(g,g') - ^ !^e'P(i-i') Exp[-{^^ + 
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peaked about qa with a width arising from quantum fluc- 
tuations (finite Mph) and from departures from equilib- 
rium (which act as an effective temperature). 

Fixing the pa requires consideration of the expo- 
nentiaUy small processes neglected above. Within the 
present formalism we find two: quantal (finite Mph) ef- 
fects which lead to tunnelling through the barrier con- 
necting the minima, and diffusion (finite A/i) effects 
which produce motion along the pseudoenergy surface. 
Diffusion may be analysed using the analogy between 
A/i and temperature mentioned above; it leads^^ to mo- 
tion between minima with rate constant Rdi f f of order 
In Riiiff « if A/i >> uJa and to a much smaller 
rate in the opposite limit. The standard equilibrium 
estimate of the rate due to tunnelling processes gives 
InRtun ~ —2^ Hi<_ a(cn-q2) ^ corrections become 
important only for A/i uja- Thus we expect roughly 
that for A// > Ua, the diffusion process dominates. 



with T ff — Au antikink. The path integral Eq|51thus becomes 
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Pait) = pait) -Rj^^ dtfdU{TZ{tf,U; {qca qca})Pa{ta) 

-n{tf ,U;{q,a ^ qcb} Pbih)) + OR^ (14) 

where the first term is the sum of the processes shown 
in panels a and b in Fig. |21 the second term is the sum 
of the processes shown in c and d, and we have denoted 
explicitly only the beginning and ending values of the 
classical component of the field. 

Differentiating Eq^]with respect to time yields 



dpajt) 
dt 



-Ta^bPait) + Tb^aPbit) 



(15) 



where to order i?(„„ the scattering rates are 



dm{t,0;{q,a^qca'}) (16) 



The integral in Eq 1161 yields a quantity Tq with di- 
mension of time; we may neglect higher order terms if 
RtunTo << 1. To evaluate Eq^|it suffices to approxi- 
mate the potential in Eq^jby the 'telegraph' form shown 
in Fig 121 the equations may then be solved by standard 
methods. In equilibrium a crucial role is played by the 
phase shift 5 = tan"! (^^S±+^}C\m ^ith V ^\{qb- qa) 

and a^R and b given by Eas I5I6I at q — qa. As shown 
by Ng^^, out of equilibrium, the behavior is described by 
complex phase shifts Sl.r with 



5l.r tan ^ 



aL.RV 



l~bV - iaR^LVsgn {pl,r - Prx) 



Note Sl + Sr^ S. We find 



(17) 



FIG. 2: Instanton processes entering equation of motion of 
density matrix. 

We now consider in more detail the interesting quan- 
tum limit A/i < Ua for large barrier height H. In this 
limit one would like to restrict attention to paths for 
which q(t) is almost always near one of the minima, with 
occasional tunnelling events in which q shifts from one 
minimum to another. The tunnelling amplitude Rtun is 
exponentially small and (for small enough A/i) is well 
approximated by its equilibrium value. However, in the 
real time path integral the smallness of Rtun arises as 
the nearly complete cancellation of a sum of many paths 
with large but oscillating amplitudes. We argue that for 
small A/i the result of performing this complicated sum 
(which we do not treat directly) may be represented by 
the 'instanton' (kink) vertices shown in Fig |21 with am- 
plitude Rtun- A kink or antikink comes with a factor i 
and a sign a — ±1 determined by the contour it is on. 
We find (see below) that the action of a single kink or 
antikink is infinite, so each kink must be followed by an 



InTZ ^ (Ceqip) - iSEneq + SCdec + SCorth)tfi (18) 

Here SEneq — (^^ — S'j^) Ap is the change in electronic 
energy arising from the imbalance in chemical potential. 
Ceqi'p) — iAEneq is the equilibrium result evaluated at 
the mean chemical potential plus the nonequilibrium en- 
ergy correction; it oscillates at frequency AE, provid- 
ing convergence of the integral in Eg 1141 when AE ^ 0. 
The imaginary parts 5" of the phase shifts give rise to a 
time decay ('decoherence') which is the manifestation, in 
the present formalism, of the decoherence introduced on 
semiphenomenological grounds by Rosch et. al.^'^. and 
ensures that out of equilibrium the integral in Eall4lcon- 
verges even when AE + AEneq = 0. We find 

SCdec = -^^ ^-ApM^ptf^) (19) 

zvr 

At short times Aptfi « 1, SCdectfi ^ ~{^ptfi)'^] 
at long times SCdectfi ~ "I^M^/il- At intermediate 
times an analytical solution is not available and (f)2 must 
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FIG. 3: Mam panel: Variation of current with chemical po- 
tential difference A/i at fixed pseudoenergy difference E1 — E2 
calculated for parameters Xqi = 4E" (resonance well above 
fermi level) and Xq2 = — E" (resonance slightly below fermi 
level) and equal coupling to left and right leads. Upper Inset: 
Spectral functions corresponding to equilibrium state (solid 
curve-only minimum " 1" occupied) and strongly nonequilib- 
rium state (dashed curve; incoherent superposition of the G's 
corresponding to the two minima). Lower Inset: Variation 
with pseudoenergy difference of density matrix element corre- 
sponding to states near pseudoenergy minimum 1, computed 
for parameters used in main panel. 



be computed perturbatively or numerically. We find 

^Gorth expresses the change in orthogonality effects 
arising because at nonvanishing A/i there is destructive 
interference between the left and right leadsii: 



SCorthtfi 

with ^"1(2;) = 
I du s\n{ux) 



2ReSLSR 2iImSLSR 



= 7e — Ci{x) + hi{x) and -02 (a;) 

ln{l-u)+u Inu] 



(20) 



crossover functions to evaluate F (Ea ll6() . In equilibrium 
and at f3j = 00 we find that only scattering from the 
higher energy to the lower energy extremum occurs (the 
'up-scattering rate' vanishes); under nonequilibrium con- 

ditions an upscattcring rate appears: of order ^-^^^ ^ 
if ^-^j^ >> 1 but of the order of the down-scattering rate 
as ^-^j^ — > 0. The change in the weight pi for one of the 
two minima as A/i/A_E is varied is shown in the lower 
inset of Fig |3| A similar result was found by ParcoUet 
and Hoolcy in a 'pseudofermion' diagrammatic study of 
the nonequilibrium Kondo problem.^'*. 

The semiclassical Greens function follows as: 



pa ((A/i) 



eo - gqa - Sfl(cj) 



(21) 



We have used the perturbative expressions for the 



The main panel of Fig |31 shows the current computed in 
the standard way by inserting the Green function given in 
Eq. 12] into Eq 28 of ReS. The inset shows the evolution 
of the spectral density A{ui) = ImGii{uj). 

The minimum condition 9$ /dq = is equivalent to the 
Hartree-Fock equations discussed in recent literature^, 
but in these works it is assumed that at each set of param- 
eter values, only one of the minima is occupied {pa = or 
1); preparation conditions are argued'^— to determine the 
state of the system, leading to multistability and switch- 
ing, in contradiction to the results shown in Fig 13 Al- 
though the bistability discussed bjs^ does not occur, the 
slow dynamics governing equilibration between the min- 
ima will lead to 'telegraph' noise in the current. 

Hamann^^ shows that interacting electron problems 
such as the Anderson model lead, after Hubbard- 
Stratonovich transformation, to a model very like that 
studied here, but with the role of the phonon field played 
by the spin part of the decoupling field. Our analysis 
carries over directly to these problems, providing a dif- 
ferent derivationi^ of the generally accepted resultsiii^ 
that nonequilibrium effects suppress the formation of the 
Kondo resonance, and that the nonequilibrium Kondo 
effect is fundamentally a weak coupling problem. 
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